# Title     : TODO
# Objective : TODO
# Created by: Administrator
# Created on: 2019/7/24

library(ggrepel)
library(optparse)
library(stringr)
library(magrittr)
library(tidyverse)
library(limma)
library(dplyr)
library(impute)

createWhenNoExist <- function(f) {
  !dir.exists(f) && dir.create(f)
}

getmode <- function(v) {
  uniqv <- unique(v)
  uniqv[which.max(tabulate(match(v, uniqv)))]
}


option_list <- list(
  make_option("--i", default = "AllMet_Raw.txt", type = "character", help = "metabolite data file"),
  make_option("--g", default = "group.txt", type = "character", help = "sample group file"),
  make_option("--config", default = "calculate_config.txt", type = "character", help = "config file")
)
opt <- parse_args(OptionParser(option_list = option_list))

args <- commandArgs(trailingOnly = F)
scriptPath <- dirname(sub("--file=", "", args[grep("--file", args)]))
source(str_c(scriptPath, "/base.R"))

configData <- read_tsv(opt$config, col_names = F) %>%
  set_colnames(c("arg", "value"))
configData

groups <- configData %>%
  filter(arg == "groups") %>%
  .$value %>%
  str_split(":") %>%
  unlist()
groups
covs <- configData %>%
  filter(arg == "covs") %>%
  .$value %>%
  str_split(":") %>%
  unlist()
covs
kinds <- configData %>%
  filter(arg == "kinds") %>%
  .$value %>%
  str_split(":") %>%
  unlist()
kinds
fillMethod <- configData %>%
  filter(arg == "fillMethod") %>%
  .$value
fillMethod
k <- configGet(configData, "knn") %>%
  as.numeric()

options(digits = 3)

mat <- read_tsv(opt$i) %>%
  rename(Metabolite = 1)
mat

colNames <- mat %>%
  select(-c("Metabolite")) %>%
  colnames()

covKindTb <- tibble(cov = covs, kind = kinds)

numberCovs <- covKindTb %>%
  filter(kind == "number") %>%
  .$cov
factorCovs <- covKindTb %>%
  filter(kind == "factor") %>%
  .$cov

sample_info_df <- read_tsv(opt$g) %>%
  arrange(Sample = factor(Sample, levels = colNames)) %>%
  mutate_at(vars(numberCovs), function(x) {
    x %>% as.numeric()
  }) %>%
  mutate_at(vars(factorCovs), function(x) {
    x %>% as.factor()
  })

sample_info_df

naSampleIds <- sample_info_df %>%
  filter_at(vars(-c("Sample", "ClassNote")), any_vars(is.na(.))) %>%
  .$Sample
naSampleIds

sample_info_df <- if (fillMethod == "delete") {
  sample_info_df %>%
    filter_at(vars(-c("Sample", "ClassNote")), function(x) {
      !is.na(x)
    })
}else if(fillMethod=="modeFill"){
  sample_info_df %>%
    mutate_at(vars(-c("Sample", "ClassNote")), function(x) {
      mode<-getmode(x)
      replace_na(x, mode)
    })
} else if(fillMethod=="fillWithMedian"){
  sample_info_df %>%
    mutate_at(vars(numberCovs), function(x) {
      median<-median(x,na.rm = T)
      replace_na(x, median)
    })
} else sample_info_df

sample_info_df

mat <-if(fillMethod == "delete"){
  mat %>%
    select(-c(naSampleIds))
}else mat

mat <-  mat %>%
  column_to_rownames("Metabolite")
mat %>%
  head()

cn <- combn(groups, 2)
for (i in 1:ncol(cn)) {
  row <- cn[, i]
  group1Name <- row[1]
  group2Name <- row[2]
  parent <- str_c(group1Name, "-", group2Name)
  createWhenNoExist(parent)

  contrast_note <- str_c("ClassNote", group1Name, "-ClassNote", group2Name)
  contrast_note

  formula <- as.formula(str_c("~0 + ClassNote + ", str_c(covs, collapse = " + ")))
  print(formula)
  design_D <- model.matrix(data = sample_info_df, formula)
  print("=in=")
  print(design_D)
  print("=out=")

  cont.matrix_D <- makeContrasts(contrast1 = contrast_note,
                                 levels = design_D)
  print(cont.matrix_D)

  fit_D <- lmFit(mat, design_D)
  fit_D <- contrasts.fit(fit_D, cont.matrix_D)
  fit_D <- eBayes(fit_D)

  Treat.adj <- topTable(fit_D, coef = contrast_note, number = nrow(mat), adjust.method = 'BH')
  tt_D <- topTable(fit_D, coef = 0)


  Treat.adj <- Treat.adj[, colnames(Treat.adj) != 'logFC']
  Treat.adj <- Treat.adj[, c('t', 'P.Value', 'adj.P.Val')]
  colnames(Treat.adj) <- c('t', "P.value", 'FDR')

  Treat.adj <- Treat.adj %>%
    rownames_to_column("Metabolite") %>%
    rowwise() %>%
    do({
      result <- as_tibble(.)
      group1SampleIds <- sample_info_df %>%
        filter(ClassNote == group1Name) %>%
        .$Sample
      group2SampleIds <- sample_info_df %>%
        filter(ClassNote == group2Name) %>%
        .$Sample
      group1Data <- mat[result$Metabolite, group1SampleIds] %>%
        unlist()
      group2Data <- mat[result$Metabolite, group2SampleIds] %>%
        unlist()
      fcList <- my.fc(group1Data, group2Data, "auto")
      result %>%
        mutate(log2FC = log(fcList$fc, 2))
    }) %>%
    ungroup()
  write_csv(Treat.adj, str_c(parent, "/res.csv"))

}














